Wang Haihua
🍈 🍉🍊 🍋 🍌
对于多元线性回归, 当 $X^{T} X$ 不是满秩矩阵时存在多个解析解, 它们都 能使得均方误差最小化, 常见的做法是引入正则化项。岭回归和 LASSO 回 归是目前最为流行的两种线性回归正则化方法。
从理论上说, 在讨论多元线性回归模型时,一般要求设计矩阵 $X$ 中的列 向量之间不存在线性关系。当 $X$ 的列向量之间有较强的线性相关性时, 即解 释变量 $x_{1}, x_{2}, \cdots, x_{m}$ 之间出现严重的多重共线性, 这时设计矩阵 $X$ 将呈病态。 在这种情况下, 用普通最小二乘法估计模型参数, 往往参数估计方差太大, 使普通最小二乘法的效果变得很不理想。为了解决这一问题, 统计学家从模 型和数据的角度考虑, 提出了很多改进方法, 下面要介绍的岭回归和 LASSO 回归是其中的两种方法。
$$\begin{array}{ccccc|ccccc} \hline \text { 年份 } & x_{1} & x_{2} & x_{3} & y & \text { 年份 } & x_{1} & x_{2} & x_{3} & y \\ \hline 1949 & 149.3 & 4.2 & 108.1 & 15.9 & 1955 & 202.1 & 2.1 & 146.0 & 22.7 \\ 1950 & 171.5 & 4.1 & 114.8 & 16.4 & 1956 & 212.4 & 5.6 & 154.1 & 26.5 \\ 1951 & 175.5 & 3.1 & 123.2 & 19.0 & 1957 & 226.1 & 5.0 & 162.3 & 28.1 \\ 1952 & 180.8 & 3.1 & 126.9 & 19.1 & 1958 & 231.9 & 5.1 & 164.3 & 27.6 \\ 1953 & 190.7 & 1.1 & 132.1 & 18.8 & 1959 & 239.0 & 0.7 & 167.6 & 26.3 \\ 1954 & 202.1 & 2.2 & 137.7 & 20.4 & & & & & \\ \hline \end{array}$$下表 是 Malinvand 于 1966 年提出的研究法国经济问题的 一组数据。所考虑的因变量为进口总额 $y$, 三个解释变量分别为: 国内总产 值 $x_{1}$, 储存量 $x_{2}$, 总消费量 $x_{3}$ (单位均为 10 亿法郎)。
对于上述问题，我们可以直接用普通的最小二乘估计建立关于三个解释变量$x_1,x_2,x_3$和$y$的回归方程为 $y=-8.6203-0.0742 x_{1}+0.5104 x_{2}+0.3116 x_{3}$ ， 并且模型的统计检验指标都相当好，但是$$x_1的系数为负，这不符合经济意义，因为法国是一个原材料进口国，当国内总产值$x_1$增加时，进口总额$y$也应该增加，所以该系数的符号应该为正。其原因就是因为三个自变量$x_1,x_2$和$x_3$之间存在多重共线性。
$x_{1}, x_{2}$ 和 $x_{3}$ 三者的相关系数矩阵 $$ R=\left[\begin{array}{ccc} 1 & -0.0329 & 0.9869 \\ -0.0329 & 1 & 0.0357 \\ 0.9869 & 0.0357 & 1 \end{array}\right] \text {, } $$ 由此可知 $x_{1}$ 与 $x_{3}$ 间的相关系数高达 $0.9869$, 这说明 $x_{1}$ 与 $x_{3}$ 基本线性相关, 若将 $x_{3}$ 看作因变量, $x_{1}$ 看作解释变量, 那么 $x_{3}$ 关于 $x_{1}$ 的一元线性回归方程为 $$ x_{3}=-4.9632+0.7297 x_{1}, $$ 这说明当 $x_{1}$ 变化时, $x_{3}$ 不可能保持一个常数, 因此对回归系数的解释就复杂 了, 不能仅从其符号上作解释, $x_{1}$ 和 $x_{3}$ 之间存在着多重共线性关系。
代码
import numpy as np; import statsmodels.api as sm
a=np.loadtxt("data/economic.txt") #加载表中x1,x2,x3,y的11行4列数据
x=a[:,:3] #提出自变量观测值矩阵
X=sm.add_constant(x) #增加第一列全部元素为1得到增广矩阵
md=sm.OLS(a[:,3],X).fit() #构建并拟合模型
b=md.params #提取所有回归系数
y=md.predict(X) #求已知自变量值的预测值
print(md.summary()) #输出模型的所有结果
print("相关系数矩阵:\n",np.corrcoef(x.T))
X1=sm.add_constant(a[:,0])
md1=sm.OLS(a[:,2],X1).fit()
print("回归系数为：",md1.params)
“多重共线性”一词由 R. Frisch 在 1934 年提出的, 它原指模型的解释变 量间存在线性关系。在实际经济问题中, 由于经济变量本身的性质, 多重共 线性是存在于计量经济学模型中的一个普遍的问题, 产生多重共线性的原 因一般有以下三种情况：
(1) 许多经济变量之间存在着相关关系, 有着共同的变化趋势, 例如 国民经济发展使国民增加了收入，随之消费、储蓄和投资出现了共同增长。 当这些变量同时进入模型后就会带来多重共线性问题。如果采用其中的两个作为解释变量, 就可能产生多重共线性问题。
(2) 在回归模型中使用滞后解释变量, 也可能产生多重共线性问题, 由于经济变量的现期值和各滞后期值往往高度相关。因此使用滞后解释变 量所形成的分布滞后模型就存在一定程度的多重共线性。
(3)样本数据也会引起多重共线性问题。根据回归模型的假设, 解释 变量是非随札变量, 由于收集的数据过窄而造成某些解释变量似乎有相同 或相反的变化趋势, 也就是说解释变量即使在总体上不存在线性关系, 其样 本也可能是线性相关的。在此意义上说多重共线性是一种样本现象。
一般地，关于多重共线性关系，我们给出如下定义。
定义1 当设计矩阵 $X$ 的列向量间具有近似的线性相关时, 即存在 不全为 0 的常数 $c_{0}, c_{1}, \cdots, c_{m}$, 使得 $c_{0}+c_{1} x_{1}+\cdots+c_{m} x_{m} \approx 0$, 称各自变量之间有多重共线性关系。
import numpy as np; import statsmodels.api as sm
a=np.loadtxt("data/economic.txt") #加载表中x1,x2,x3,y的11行4列数据
x=a[:,:3] #提出自变量观测值矩阵
X=sm.add_constant(x) #增加第一列全部元素为1得到增广矩阵
md=sm.OLS(a[:,3],X).fit() #构建并拟合模型
b=md.params #提取所有回归系数
y=md.predict(X) #求已知自变量值的预测值
print(md.summary()) #输出模型的所有结果
print("相关系数矩阵:\n",np.corrcoef(x.T))
X1=sm.add_constant(a[:,0])
md1=sm.OLS(a[:,2],X1).fit()
print("回归系数为：",md1.params)
OLS Regression Results
==============================================================================
Dep. Variable: y R-squared: 0.996
Model: OLS Adj. R-squared: 0.994
Method: Least Squares F-statistic: 542.0
Date: Sun, 22 May 2022 Prob (F-statistic): 1.20e-08
Time: 14:28:48 Log-Likelihood: -1.7480
No. Observations: 11 AIC: 11.50
Df Residuals: 7 BIC: 13.09
Df Model: 3
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t| [0.025 0.975]
------------------------------------------------------------------------------
const -8.6203 0.897 -9.611 0.000 -10.741 -6.499
x1 -0.0742 0.028 -2.691 0.031 -0.139 -0.009
x2 0.5104 0.075 6.781 0.000 0.332 0.688
x3 0.3116 0.037 8.357 0.000 0.223 0.400
==============================================================================
Omnibus: 5.258 Durbin-Watson: 2.400
Prob(Omnibus): 0.072 Jarque-Bera (JB): 2.250
Skew: 1.080 Prob(JB): 0.325
Kurtosis: 3.495 Cond. No. 2.05e+03
==============================================================================
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
[2] The condition number is large, 2.05e+03. This might indicate that there are
strong multicollinearity or other numerical problems.
相关系数矩阵:
[[ 1. -0.03291436 0.98690551]
[-0.03291436 1. 0.03567322]
[ 0.98690551 0.03567322 1. ]]
回归系数为： [-4.96322784 0.72966696]
D:\software_install\anaconda\lib\site-packages\scipy\stats\stats.py:1541: UserWarning: kurtosistest only valid for n>=20 ... continuing anyway, n=11
warnings.warn("kurtosistest only valid for n>=20 ... continuing "
a
array([[149.3, 4.2, 108.1, 15.9],
[171.5, 4.1, 114.8, 16.4],
[175.5, 3.1, 123.2, 19. ],
[180.8, 3.1, 126.9, 19.1],
[190.7, 1.1, 132.1, 18.8],
[202.1, 2.2, 137.7, 20.4],
[202.1, 2.1, 146. , 22.7],
[212.4, 5.6, 154.1, 26.5],
[226.1, 5. , 162.3, 28.1],
[231.9, 5.1, 164.3, 27.6],
[239. , 0.7, 167.6, 26.3]])